


Institutional Archive of the Naval Postgraduate School 





Calhoun: The NPS Institutional Archive 
DSpace Repository 


Theses and Dissertations l. Thesis and Dissertation Collection, all items 


1972-06 


Kalman filtering techniques applied to 
airborne direction-finding and emitter location 


Colburn, L. Laddie. 


http://hdl.handle.net/10945/16229 


This publication is a work of the U.S. Government as defined in Title 17, United 
States Code, Section 101. Copyright protection is not available for this work in the 
United States. 


Downloaded from NPS Archive: Calhoun 


Calhoun is the Naval Postgraduate School's public access digital repository for 
| г D U DLEY research materials and institutional publications created by the NPS community. 
FW узу, Calhoun is named for Professor of Mathematics Guy K. Calhoun, NPS's first 


ІШІ KNOX appointed — and published — scholarly author. 

OM LIBRARY Dudley Knox Library / Naval Postgraduate School 

411 Dyer Road / 1 University Circle 
Monterey, California USA 93943 





http://www.nps.edu/library 


KALMAN FILTERING TECHNIQUES APPLIED TO 
AIRBORNE DIRECTION-FINDING AND 
EMITTER LOCATION 





L. Laddie Coburn 


a mein 











MES 
τσεχ 
«гг E 
ЖС 
ᾱ--- <x... 
c3 
exo") 
ο) 
c s 
acto NOT е 
Bp 5 ТАУ 
> 
a 
τ 
с a 
p a > = 





Іі а - “€ i ο 
| | S í ° 


x KALMAN FILTERING TECHNIQUES © 
APPLIED TO 


AIR BORNE DIRECTION-FINDING 
AND EMITTER LOCATION 


7, hess s 
UNS by 


L. Laddie Coburn 





Thesis Advisor: H. A. Titus 
June 1972 


Approved for public release; distribution unlimited 





De 





ADUATE SCHOOL 


que 
шакта 


Monterey, Galifornia 










KALMAN FILTERING TECHNIQUES 
APPLIED TO 
AIR BORNE DIRECTION-FINDING 
AND EMITTER LOCATION 


by 


L. Laddie Coburn 






Thesis Advisor: A ИЕ 
| June 1972 


Approved for public release; distribution unlimited 





Kalman Filtering Techniques 
Applied to 
Airborne Direction-Finding 
and Emitter Location 


by 


L. Laddie Coburn 
Lieutenant, United States Navy 
B.S., United States Naval Academy, 1965 


Submitted in partial fulfillment of the 
requirements for the degree of 


AERONAUTICAL ENGINEER 


from the 


NAVAL POSTGRADUATE SCHOOL 
June 1972 





ἐστι C | 


The objective of this analytical study was to develop an optimal 
emitter location algorithm using Kalman filtering techniques to filter 
emitter bearing angies-of-arrival information in an airborne, multi- 
emitter environment. Since emitter and aircraft navigation data were 
sampled at discrete time intervals, discrete sequential estimation 
techniques were utilized within a software computer program which 
initially sorted data, filtered emitter bearing angles-of-arrival to 
obtain optimal estimates, and computed emitter positions in latitude/ 
longitude coordinates using various emitter location algorithms. An 
error analysis was performed onthe emitter location algorithms to 
determine relative accuracy oí the different techniques. Extended 
Kalman filtering was also considered as an alternate approach to the 


emitter location algorithm. 
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1. INTRODUCTION 


One goal of airborne direction-finding (DF) is to obtain the 
location of an electromagnetic radiation source in latitude/longitude 
coordinates on the earth's surface. Over the years many emitter 
location techniques have been devised, and much research has been 
devoted to optimizing bearing angle measurements in angle-only 
direction-finding and navigation systems. Much of the effort has been 
directed toward physical components such as the directional receiving 
antenna and the improvement of its characteristics. In airborne DF 
systems current typical bearing angle-of-arrival (AOA) accuracy is 
found to be about + 2° at microwave frequencies. Despite such 
accuracy a multiple bearing fix will yield considerable emitter position 
ambiguity at any significant range between the emitter and the airborne 
DF system. 

To improve this situation a statistical optimal estimate of each 
received emitter bearing angle was computed using sequential estima- 
tion techniques, commonly called Kalman filtering. Additionally, an 
adaptive gating technique was utilized which selectively filtered out 
bearing lines which did not closely associate with other bearing lines 
already filtered and correlated to a distinct emitter location. Smooth- 
ing techniques were utilized to improve the unfiltered initial bearing 


angle, and a smoothed initial bearing angle and filtered final bearing 





angle were used in a spherical earth triangulation solution to estimate 
the emitter position. The entire sorting, filtering, smoothing, gating, 
and triangulation procedure was organized within a software computer 
program capable of automatically estimating multiple emitter locations 


given noisy emitter bearing angle data and aircraft navigation data. 





Ii, ANALYSIS OF THE PROBLEM 


A. APPLICATION 

The objective of this analytical study was to develop an optimal 
emitter location algorithm which would have абым да to military 
aircraft flights in an operational environment. The computer program 
developed could be utilized in an airborne digital computer system to 
give real-time analysis of emitter parameters in a multiple emitter 
environment, or it could be used in a ground-based computer system 
to give post-flight analysis of large quantities of compiled data. 

Since all data sampling and processing was done at discrete time 
intervals, discrete sequential estimation techniques were utilized. 

The real facility in using sequential estimation techniques is that each 
subsequent computation is based only on the new observation and the 
last previously calculated estimate, and these recursive equation com- 
putations require a minimum of computer storage to filter a large 
quantity of data. This makes such an emitter location algorithm ideal 
for airborne computer application [1]. 

To allow generality of application, the software system was 
developed to process multiple emitter locations, to sort and filter 
large quantities of emitter data, to filter data at non-uniform sampling 
intervals, and to estimate emitter locations on any region of the earth's 


surface. It was initially assumed that no errors existed in the aircraft 





navigation data; hence these position data together with two or more 
filtered emitter bearing angles could be used to compute the emitter 
location (see Figure 1). If a deterministic aircraft navigation error 
was known to exist, it could be compensated for in the emitter location 
algorithm. If a random error existed in the aircraft navigation data, 
it could be minimized by using a non-zero Q matrix to account for 
random excitation noise (see Section III, Part B, page 17) or by Kal- 
man filtering the aircraft navigation data as well as the emitter bear- 
ing data. Kalman filtering techniques have also been developed to 


optimize navigation system accuracy [2]. 
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FIGURE 1 


GEOMETRY OF EMITTER LOCATION 





В. KALMAN FILTER UTILIZATION 

A signal X, measured in the presence of random white Gaussian 
noise V, is filtered in an attempt to recover the original signal. The 
filtered estimate of the state is denoted by X. From Kalman filter 
theory [3] for discrete systems using first-order difference equations, 


the kth noisy observation Z(k) of the state is given by 


ντ У Ke > pee eer yt os (1) 
where H(k) is the observation matrix and V(k) is the measurement 
noise. Since a discrete or sampled system is being considered, k 
denotes the kth time sample. 

The state of a general, discrete system with no deterministic 
forcing function is given by 

X(k+1) = O(k+1,k)X(k) + IT W(k) κου ο... (9) 
where Ф (к+1, k) is the state transition matrix, T is a distribution 
matrix related to the random forcing function, and W(k) is a random 
forcing function to account for random excitation noise. It was 


assumed that the noise sequences have zero mean and second-order 


statistics described by 


E[ vaov | = R(k)8(kj) (3) 
Ге[уо ц? | T? = Qtk) Sj) (4) 
and E| vaow? ] = 0 for all ik, j | | (5) 
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0 kzj 
where O (kj) = is the Kronecker delta. 
1 k=j 


The Kalman filter recursion equations [1] are summarized 
aS 
below where X(k/j) denotes the estimate of the state X(k) based upon j 


measurement observations Z(1), Z(2),...., Z(j). 


G (k) = р(к/к-1)н(®)Т | H(k)P (k/k-1)H x)" +в] (6) 
БЕН К CI PE Op 1) (7) 
р(к+1/к) = Ф(к+1, к)р(к/к) Ф(к+1, к)Т + ©(к) (8) 
ZN ZN IN 

Х(к/к) - Хк/к-1) е о) 2Хк) - н(Ә/Х(к/к-1)) (9) 
Sk/k-1) = D(k, k-1)X(k-1/k-1) (10) 


The Kalman filter gains are represented by the matrix G(k), and 
P(k/j) represents the error covariance matrix of the various estimates. 
For the problem being considered the states are bearing angle- 
of-arrival and bearing rate, and it is desired to filter the noisy 
observed bearing angles to obtain optimal estimates of these states. 
The system dynamics were approximated by a 1/52 plant which results 
in a state transition matrix given by 
ІР ТЕТІ) 
Ф(к-1,к) = 
0 1 
Since only the bearing angle of the state was observed, R(k) and W(k) 
become scalar variance terms and the bracketed terms in equations 
(6) and (9) become scalars as well. This considerably simplified 
reducing these matrix equations to scalar form (see Appendix A); 


however, allowing a non-uniform sampling interval required that the 
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term T(k+1) remain a variable in each of the applicable recursion equa- 
tions. This required that error covariance terms and Kalman filter 
gains be computed on line since they were a function of the sampling 
interval; therefore, these terms could not be precomputed and stored. 
Several assumptions are implicit here. Modeling the system 
dynamics with a 1/s? transfer function gives a linear approximation 
that may not adequately represent the actual system equations, which 
are non-linear [4]. It may then become necessary to augment the 
dynamic model with an additional difference equation which includes 
the second derivative of the measured bearing angle, or extended 
Kalman filtering techniques may be employed to more accurately 
approach the non-linearity of the dynamic system. Extended Kalman 
filtering techniques are discussed in Section III, Part F, of this paper. 
The statistical nature of Kalman filtering techniques requires 
that the statistics of the measurement noise V(k) and the excitation 
noise W(k) be completely specified, or if unknown they must be 
assumed. If no deterministic errors are known to be coloring these 
noise sources, the assumption of white noise may be a good one al- 
though it may not exactly model the actual noise bearing error. Kal- 
man filtering assumes that the noise is independent from one sampling 
interval to the next [5]. This may not be an accurate assumption since 
successive bearing angle observations from a single emitter might 
have similar error statistics. Any deterministic bias in aircraft 
antenna bearing angle-of-arrival measurement would further negate 
this assumption. This difficulty was reo by Reeves [6], who recom- 


mended that emitter bearing angle-of-arrival measurement noise be 
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considered to have both random and correlated components. This 
would require augmenting the state vector equation with an additional 


variable and further complicate the Kalman filtering process. 


C. SMOOTHING TECHNIQUES 

Since Kalman filtering techniques sequentially filter the second 
and subsequent bearing angle observations associated with a single 
emitter, the first bearing angle associated with a particular emitter 
remains unfiltered and no optimal estimate exists. Thus a technique 
of filtering backward in time or smoothing was used to improve the 
estimate of the initial bearing angle. The smoothing equations are 
similar to Kalman filter equations in that they are recursive and 
functions of similar statistical parameters, but they update or smocth 
previous data based on more recent data that has been optimally 
estimated. 

A general technique proposed by Rauch [7] gives the smoothed 


estimate of the initial observation after k time intervals as 


IN IN IN 
&a hio » Sa k-) + ра люс [200 - но /s-0] (11) 
where D(1/k) » D(1/k-1)P(k-1/k-1) (x, k-1) T P(c/k-1): (12) 


For a non-uniform sampling interval D(1/k) becomes a function of 
T(k) and must be solved sequentially before the smoothed initial 
estimate 1 /k) can be computed. For an analysis of these equations 
reduced to scalar form, see Appendix B. These scalar smoothing 
equations were included in the software computer program and yielded 


a smoothed estimate of the initial bearing angle and bearing rate. 
ES 





H. COMPUTATIONAL PROCEDURE 


A. INITIAL DATA SORT 

The problem was posed to accommodate a realistic airborne 
environment with multiple emitters at various RF frequencies, emitter 
bearing angle-of-arrival errors, and a constant velocity moving air- 
mait. Aircraft navigation and emitter target data was sampled and 
recorded at discrete but time-varying intervals. Since the number of 
emitters was unknown and large quantitites of signal data were recorded 
at various frequencies, pulse repetition frequencies (PRF), and pulse 
widths (PW), it was necessary that the data be initially sorted by fre- 
quency and PRF to initially estimate the number of distinct emitter 
targets. This was accomplished with sufficiently small frequency and 
PRF gates to separate and associate emitter parameters with a distinct 
emitter target. Additionally, it was found that multiples of PRF may 
exist which associate with a single emitter, so a check for PRF multi- 
ples had to be included within the PRF Sorting subroutine. One set of 
data processed by this program included PRF multiples given on page 99. 

Multiple emitter 5 of the same frequency and PRF could still 
result in an ambiguous solution; therefore, each bearing angle initially 
associated with a distinct emitter by the initial data sort was sequen- 
- tially filtered and gated to determine if it COS with this emitter. 


Due to the large quantity of data to be processed, any bearing angle not 
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correlated to a distinct emitter after being filtered and gated was dis- 
carded from the data bank. Further investigation and preater computer 
Storage Capability could allow such discarded bearing angles to be re- | 
associated with a different emitter; however, this process was not 


developed within this analysis. 


B. KALMAN FILTER INITIALIZATION 

Initialization of the Kalman filter equations requires a knowledge 
of the system dynamics, the statistical properties of the Kalman filter 
parameters, and some information about the initial state of the system. 
Since it was necessary to associate each sorted data sample (emitter 
bearing angle and associated emitter ІНІ with a particular emitter. 
the notation JSET (I, J) was devised, where I is the emitter target 
number with which the data sample is associated and J is the jth sample 
of sequential data associated with that emitter. JSET(I, J) gives the 
Sequential sample number for this data sample when listed together 
with all other data. See page 63 for a listing of initial JSET data. To 
initialize the Kalman filter equations, KI was used as the initial JSET 
value; i.e., JSET(I, 1) = KI for the ith emitter target. 

To filter emitter bearing angles-of-arrival, the noisy observation 
of the state from equation (1) is given by 


Ө(к) - Н(к) Ө(к) + V(k) (13) 
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O (k) 


where ACI = vector of noisy measured observations, 
6 (19) 
@(k) 
@(k)=| , = state vector of exact emitter bearing angle 
@(k) 


and bearing rate, and H(k) = [1 0 | is the observation matrix for 
measurement of only noisy bearing angle @. Since the bearing rate 

0 (k) was not measured, it must be estimated to initialize the sequen- 
tial estimation process. This bearing rate depends on the speed and 
heading of the aircraft with respect to the emitter, the relative bearing 
angle-of-arrival, and the range r to the emitter. Since the range to an 


unknown emitter is unknown, it must be estimated. For this program 


r = 150 nautical miles was assumed and 


Q(KI) = Z sin (BRNG) x PRAD 


Ж рс (14) 
where BRNG = relative bearing angle-of-arrival 

PIRAD = 57.29578 degrees/radian 

v = aircraft velocity in knots. 
If an estimated value of range to the emitter is known and different 
from 150 nautical miles, then this new value should be used in 
equation (14) and substituted into the computer program. 

Since an initial optimal estimate was not available to the filtering 

process described by equations (9) and (10), the first bearing angle 


optimal estimate THTD(KI) was assumed equal to the noisy measured 


bearing angle THETAD(KI) and the first bearing rate estimate TDTD(KI) 
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was assumed equivalent to (KD. The smoothing process was simi- 
larly initialized. 

The initial uncertainty of emitter bearing angle and angle rate on 
filter initialization was accounted for in the initial values of the error 
covariance matrix P(1/0). Since stationary emitter locations were 
considered, the random forcing function W(k) of the dynamic state 
equation (2) could be set to zero if exact knowledge of the observer's 
aircraft position were known. However, since the aircrait position 
may have both random and deterministic errors associated with its 
position measurement and because of the simplified Ф matrix, a non- 
zero Q matrix was utilized. For W(k) having a zero mean and 


Ewa | - 1.0 degrees? from equation (4) 











T^() | т4(к) Т3(к) 
2 2 4 2 
Q(k) - ales ao | : ES 
T(k) EDS πιο. 
2 


The value of R(k), the scalar variance of the measurement noise, 
was assumed to be constant and estimated according to the angle meas- 
urement accuracy of the DI" system being considered. The initial D 
matrix for the sequential smoothing process was initialized as the 
identity matrix (see Appendix B). 

When utilizing this program for emitter location with different 
known measurement and covariance statistics, these initialization 
parameters must be changed in the computer program to ensure optimal 


solution locations. The only other term which must be modified by the 
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program user is the value of NUM, which is the total number of data 
samples to be sorted and filtered by this program. Since the data for 
each emitter must be filtered Separately, the filter process occurs 

once for each distinct emitter target as determined by the initial data 
sort. Once an optimal estimate of emitter bearing angle and bearing 


rate have been initialized, the filter process is ready to commence. 


om KALMAN FILTER PROCESS 
The Kalman filter equations utilized to filter the states, emitter 
bearing angle 8 and bearing rate Q , are those noted on page 11 shown as 
- ^N 
equations (6) through (10) with X replaced by dr Since these are vector 
equations and tedious to handle within a computer program, their 
scalar counterparts were derived (see Appendix A) and utilized within 


the program. The optimal estimation equation is similar to equations 


(9) and (10) and is given by 


N ÍN 
B(k/k) JL MINOR (k-1/k-1) G1(k) 
Ж ш Αν == E(k) (16) 
B(k/k) oO. 1 B(k-1/k-1) G2(k) 
AN AN 
where E(k) = 0) - 8(k-1/k-1) - T(k)8(k-1/k-1) (17) 


results in a scalar term since only one of the state variables was 


measured. 

Due to the non-uniform sampling interval, the Kalman filter gains 
Gi(k) and G2(k) may vary and not approach a steady state value at a 
uniform rate as they would with a constant sampling interval. This is 


Shown graphically in Figure 2, which depicts the transient Kalman - 
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KALMAN FILTER GAIN SCHEDULE 
FOR A NON-UNIFORM SAMPLING INTERVAL 
KALMAN FILTER PROCESS 


1.0 t 


0.8 


© 
° 
Ον 


т KALMAN FILTER GAIN 
+ t + + 
+ + Ba + 
е . + 


KALMAN FILTER GAIN G (K) 
© 
T 


0.2 


0.0 - 
9 100 200 300 &00 500 600 
ELAPSED TIME t(k) IN SECONDS 


NOTES: 


1. Points plotted above represent ‘discrete values of Kalman filter gain 
G1(K) associated with filtering the bearing angle of the state vector. 


2. Since the prediction covariance matrix (P(k+1/k) given by equation (8) 
is a function of Ф(к+1, k), which contains the .sampling interval T(k+1), 

it can be seen that the Kalman filter gain G(k) will also vary with T(k+1) 
as in equation (6).. As the sampling interval increases, the gain will also 
increase. | 
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FIGURE 3 


KALMAN FILTER GAIN SCHEDULE 
FOR A UNIFORM SAMPLING INTERVAL 


1.04 + KALMAN FILTER PROCESS 
+ 
0.8 
+ 

š 
6.6 + 
© 
< + 
5 + 
a ar i KAIMAN FILTER GAIN 
Ho ++++++ +++ +++ + +++ +++ + 
: 
3 

0.2 

0.9 

9 60 129 180 240 300 360 
ELAPSED TIME t(k) IN SECONDS 

NOTES: 


1. Points plotted above represent discrete values of Kalman filter gain G1(K) 
associated with filtering the bearing angle of the state vector. 


4. After ten sampling intervals the gain has nearly reached a steady state 
value of 0. 417. This asymptote is determined by the non-zero value of the 
Q matrix. 


3. If system random excitation noise (i.e., the Q matrix) were assumed 
to be zero, then the Kalman filter gain would approach zero upon reaching 
Steady state. 
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filter gain schedule of a non-uniform sampling interval. Figure 3 
shows the Kalman filter gain schedule of a uniform, twelve-second 
sampling interval computer simulation, which was developed to allow 
an error analysis to be made on the Kalman filter prograrn. 

The adaptive gate was chosen to be а function of the prediction 
variance term P11(k/k-1), which is a "DS function, and the 
measurement noise variance R(k), which is constant. The adaptive 
gate was set to selectively compare the filtered bearing angle with 
previously correlated bearing angles. If a filtered bearing does not 
fall within the adaptive gate, it 1S discarded and the sampling interval 
increases to the "ms previous correlated bearing. This increases the 
Kalman filter gain and opens the adaptive gate to a wider allowable 
value until a subsequent bearing angle is found which correlates to this 
emitter, then causing the adaptive gate to reduce in size. The adaptive 


gate was chosen to be 





GATE(K) = C+ «4 P1100 * ROO | (18) 
where C is a constant scale factor. Since the prediction covariance 
P(k/k-1) may be defined as 

р(к/к-1) = E Í| 80) - Su] [екю - @кук-1)| }, аз) 
which is a function of the mean square error between the true value of 
the bearing angle and its cptimal estimate, it can be seen that as the 


estimate error increases the adaptive gate wil! increase and as tre 


optimal estimate approaches the true value the gate size will decrease 
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as is desired. The constant term R(k) prevents the gate from getting 
so narrow that no emitter bearing angles could correlate to a distinct 
emitter, and the variable term P11(k) allows the gate size to "adapt" 


to the transient filter errors as was noted above. 


D. INITIAL BEARING ANGLE SMOOTHING PROCESS 

In order to obtain a smoothed estimate of the first bearing angle 
associated with each emitter target, it was necessary to solve equation 
(11) sequentially within the Kalman filter process since the gain term 
was recomputed after each sampling interval. Тһе smoothed initial 
bearing angle estimation equation used in the software program is 


similar to equation (11) and is given by 


DEN ZN 

B(1/k) B(1/k-1) D110/1) DI2CpO 9) C 

ж AS + E(k), (20) 
B(1/k) 0(1/к-1) D211 /kyY D22 E 


where E(k) is defined by equation (17) and the scalar terms in the D 
matrix are as derived in Appendix B. To compute these terms, both 
@(k,k-1)7) and P(k/k-1)7! had to be determined. The inverse of the 
state transition matrix is obviously given by 
1] -T(k) 
Ф(к,к-1)71 = . (21) 
0 ] 
Should the error covariance matrix become singular because of numer- 
ical inaccuracy in computation, no inverse will exist. In this case 


corrective action must be taken before attempting to compute the 


inverse. Several references have discussed this problem, suggested 
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solutions, and Schmidt [8] gives a survey of current techniques in this 
area. For this program, the IBM 360 subroutine MINV was utilized to 
compute the error covariance matrix inverse with good numerical 


accuracy. 


E. EMITTER LOCATION ALGORITHMS 

several techniques were investigated to obtain a computationally 
simple yet accurate algorithm to calculate the emitter locations. Since 
more than two bearing lines may result in an ambiguous emitter loca- 
tion, it was most logical to use only the smoothed initial bearing angle 
А ^N 
0; and the filtered final bearing angle 0; to compute the location. Any 
such location procedure would give most accurate results if the two 
bearings intersected at right angles; however, infrequent data collec- 
tion, short-leg flight tracks, or long-range distances from the emitter 
generally prevent such choice of initial and final bearing angles. 

1. Plane Triangulation Solution 

The basic "flat earth” triangulation solution was presented 

by Kayton and Fried [9] and discussed in some detail in Refs. 6 and 10. 


Figure 1 depicts the geometry of the triangulation method of emitter 


location solution, and the solution equations are given by 


X x 
ж PE 
tan 1 = уу (22) 
Т 1 
X X 
^A Е 
апа tan Bf = en 
Ya = Ye (23) 


where the emitter location coordinates (Xm, Y y) are unknown and must 
be determined. 


23 





If latitude (L)/longitude (X) coordinates are substituted for 


the X/Y axes in the flat earth model, the following equations result: 


MN (Ат -A,)cosL,, 


tan, = (24) 
Lp = L; 
- L 
| (A, Af)cos T 
tan, = (25) 
Lip ту L+ 


These equations gave a more accurate solution but were more complex 
to program Since one of the unknowns (Т) occurs both as an explicit 
function and as a coSine function. This difficulty was resolved by 
solving for an intermediate value of emitter latitude (TILA) based on an 
aircraft midlatitude cosine function (WAV) which was known and then by 
recomputing emitter latitude (TLA) based on the cosine function of the 


intermediate emitter latitude. The resulting equations are given by 


L. + Le 


WAV = Bur | (26) 


ж ж 
(А, - Aj)cos(WAV) * L.tang. - L,tan®, 
RENS er Se. шыс к= (27) 


Ж Ж 
tanB; - tandr 


^ Αν 
(Af -A,)cos(TILA) + L,tanB. - L tan 0, 


TLA sS Beer c c EE (28) 
tan; - бап, 
u ЖА 
TLO = A; - (Lip - І.,Жал0,/сов(ТП.А) (29) 


Т 


l! 


where TLA = emitter target latitude L, 


TLO = emitter target longitude Ay 
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L. = initial аїгеоган осол ае 


À; = initial aircraft position longitude 
Le = final aircraft position latitude 
Af = final aircraft position longitude. 


These equations are derived in Appendix С and are equivalent to the 
"round world" model developed on pages 41-43 of Ref. 6. 
2. Spherical Triangulation Solution 
Two approaches to a ο μαι earth solution were taken. 
The first involved the derivation of a spherical triangulation solution 
with two equations resulting from two spherical triangles as shown in 
Fig. 4. These equations are derived in Appendix D and are given by 


7а tan | us - d ,)cos Lo 


tan; = (30) 
sin(Lrp x L.) 
tan (gu o - A.)cosL 
сап, = tan [Ap - Apeost., | (31) 


sin(L,, - L.) 

These equations represent exact solutions for emitter loca- 
tion on a spherical earth of constant radius. It is known that the earth 
is an oblate spheroid so small errors may be introduced at large dis- 
tances from the emitter with such a spherical earth solution. However, 
according to Kayton and Fried [9], a constant radius of curvature can be 
used everywhere on the earth with an error of less than 0.2 percent of 
distance. 

If the emitter location problem is concerned only with micro- 


wave frequency emitters, then relatively short distances from the 
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aircraft to the emitter will exist since microwave elecir oma mene 
(EM) energy propagation approximates a straight line path (line-of- 
sight) on the earth's surface. Reintjes and Coate [11] show that the 
maximum distance for direct EM wave propagation to a radar horizon 


js given by 


d -«2h, * 42h, , | (32) 
where d= radar horizon in statute miles 

ha = aircraft altitude in feet 

ρ = emitter altitude in feet 


assuming a 4/3 effective radius of the earth for standard EM energy 
refraction characteristics. If h, is unknown or assumed to be zero, 
then aircraft altitudes (h,) up to 80,000 feet should not yield radar 
horizons greater than 400 statute miles (350 nautical miles). At this 
distance latitude/longitude differences between the aircraft and emitter 
are small (less than 6°), and it can be seen that the emitter location 
algorithm developed from equations (24) and (25) is nearly equivalent 
to the exact spherical solution from equations (30) and (31) since the 
tangent and sine functions of small angles are approximately equal to 
the radian measure of the angles themselves. Thus for ranges up to 
360 nautical miles from the emitter, the planar solution will be less 
than 0.5 percent in error from the spherical solution. For ranges of 
less than 150 nautical miles the error may be neglected and the plane 
triangulation solution may be considered sufficiently accurate for 


emitter location estimation. 
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If greater accuracy is desired, the 5рһег1са н тапа ы ын 
equations (30) and (31) may be solved by an iterative procedure as 
suggested in Appendix E. It was found to be tedious to attempt a 
closed-form solution. 

3.  Earth-Centered Coordinate System Solution 

Another approach to the spherical kasd solution was con- 

sidered using tangent plane projections and an earth-centered coordi- 


nate system (x, Χρ» Χα) where the x, axis is the polar axis and the xj, 


3 
X9 plane is the equatorial plane of the earth. An observation plane was 
defined for each emitter-bearing line, which passed through the center 
of the earth and inscribed a great circle arc on the earth's surface 
between the aircraft and emitter positions. The procedure for obtain- 
ing the emitter location required solving the equations of the SR 
tion planes, two at a time, with a constraint equation that the emitter 


be located on the earth's surface. Assuming a constant earth radius 


and letting X4 equal a constant, the emitter position may be computed 


from the equations of the plane (хү, хо, с) where x, = constant as shown 
below. 
, "EL 
а11Х1 + ajo9X9 = 8190 (33) 
, ^ 
ao 1*1 ЕЕ 299х9 = -а93с (34) 


These equations were then transformed back to the original 


coordinate set (%1,%9:%,) by the equation 


x, = AX, | ντ ο ο. (55) 





where А = Re/ (x ° i хо“ + ο (36) 
and Re = average radius of the earth. 

With proper coordinate rotation, it is possible to have the projection 
plane (xg = ο) tangent to the earth's surface and centered at the aircraft 


observation position. 


Ек, EXTENDED KALMAN FILTERING 

It may be desirable to directly filter and update the estimates of 
emitter position (Lm, Ат) ав described by Reeves [6]. If only bearing 
angle-of-arrival information is available as observable data, however, 
a non-linear transformation must be utilized to transform observable 
bearing angles into filtered position estimates. 

This was accomplished by extended Kalman filtering techniques 
utilizing the concepts of small perturbation theory and a Taylor series 
expansion about an initial point. Equation (1) then becomes 

ζω = м|ха) + vw { (87) 
where Z(k) represents the observable bearing angles, X(k) is the 
emitter position vector, and N represents the non-linear transforma- 


tion given by 


(A, - A)cosL 


IN -1 T 


8 = tan (38) 
Lip zu 
If the position error of the state vector X(k) is given by 
X(k) = х(к)-%(), (39) 


then the true position X(k) may be expanded about the most recent 


^^ 
optimal estimate X(k) in a Taylor series expansion with higher order 
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terms neglected as shown below: 


N [xa] = x [S09 | + MIX +22... (40) 
IN 30 ^ 

where WIR) геа = Κα | (41) 
Se E 


for the emitter location algorithm given by equations (24) and (25). 
The Kalman filter recursion equations (6) through (8) may then be 


rewritten to account for this non-linear observation matrix and are 


given by 
бк) = Р(к/к-1)м(к)Т | М(К)Р(К/к-1)М(К)Т + R0o) | (42) 
P(k/k) = P(k/k-1) - G(k)M(k)P(k/k-1) (43) 
р(к+1/к) = Ф(к+1, к)р(к/к) Ф(к+1, к) + Q(K) . (44) 


Since the non-linear filter is processing estimates of a stationary 
emitter position, the Ф matrix becomes the identity matrix and the 
equations are considerably simplified and reduced to scalar form 
(see Appendix F). 

The greatest difficulty in this approach was found in initialization 
of the extended Kalman filter process. Errors or inaccuracies in 
initialization may cause the recursive computations to diverge, and it 
was therefore decided to commence extended Kalman filter processing 
only after linear filtering and smoothing of the bearing angle observa- 
tions had been completed and an initial optimal estimate of emitter 
position obtained. This would allow the СІП Kalman filter process 
to be used as a final correction procedure to furtter improve tt#emm и 


optimal estimates of emitter position. 
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IV. PRESENTATION OF RESULTS 


A. MONTE CARLO SIMULATION ANALYSIS 
Since this Kalman filter program was developed to process actual 
aircraft flight data, a Monte Carlo simulation — utilized to accom- 
plish program validity checks and error analysis studies. Figure 5 
depicts the geometry of the simulation with an aircraft flying due south 
at 600 knots off the coast of southern California. Two VORTAC stations 
were chosen as emitters with accurately known positions given A 
EMITTER #1 Oceanside VORTAC 392240559 N 
Channel 100 -117.41694° W 
Fremguency' 119" mtijz 
EMITTER #2 nan Diego VORTAC EO DEN 
Channel 125 -117.22444? W. 
Frequency 1212 mHz 
Bearing angle-of-arrival information was assumed to arrive every Six 
seconds alternately from each target giving a twelve-Second uniform 
sampling interval for simplicity. Known true bearing angles were 
computed numerically in the absence of meaSurement noise and allowed 
— location within twenty feet of the known locations given above. 
In the presence of computer-generated random measurement 
noise of one degree variance, a large number of simulation runs were 
executed to obtain a statistically valid error analysis. The results of 


this analysis are given in Appendix G and summarized below for the 


two simulation emitters. A range error was computed from each 
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FIGURE 5 
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position estimation latitude and longitude error, and all errors were 


statistically averaged for the total number of runs. 


MONTE CARLO SIMULATION ERROR ANALYSIS 


average error in nautical miles 


normalized 

average average average closest estimated 
latitude longitude range point of position 

emitter error error _error approach range error 
1 0.44852 ο δυο 0 55692 54.5 0% 985% 
2 059 01 0229108 0950798 64.7 0. 785% 


These results show that the basic plane triangulation solution emitter 
location algorithm is capable of estimating position within one percent 


of range at emitter distances of less than 100 nautical miles. 


B. AIRCRAFT FLIGHT DATA ANALYSIS 

The Kalman filter program was also utilized io process andil! r 
aircraft flight data of the form discussed in Section III, Part A, of this 
paper. This data processing tested maximum utility of the computer 
program since certain data samples would not correlate with predicted 
emitter locations and had to be selectively gated out. Figures 6 through 
10 show graphs of noisy unfiltered bearing angle-of-arrival information 
initially associated with a distinct emitter target. The dotted lines show 
unfiltered bearing angles and the solid lines show the smoothed initial 
and filtered final bearing angle. It can be seen that Figures 7, 9, and 


10 have erroneous bearing angles which must be excluded from the 
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Kalman filtering process by the adaptive gating technique. This is 
shown to occur in the results of the computer output on page 64. Also 
comparing the final JSET data (page 65) with the initial JSET data 
(page 63) shows which data samples were not correlated to a distinct 
emitter target and were discarded. 

The emitter location estimation accuracy obtained from the 
results of processing and filtering flight data was not as good as that 
obtained from the Monte Carlo simulation; however, it was estimated 
that the standard deviation of bearing angle measurement noise was 
approximately 12° and that some aircraft navigation error may have 
existed as well. Therefore, typical emitter location estimation within 
3 percent on an average range from aircraft to emitter target of 100 


nautical miles was considered reasonable. 
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FIGURE 6 
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FIGURE 8 
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FIGURE 9 
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FIGURE 10 
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V. CONCLUSIONS 


It was found that Kalman filtering techniques have useful applica- 
tion in the field of airborne direction-finding and emitter location pro- 
cedures. Many airborne directional antenna receiving systems have 
+10° typical bearing s: which yield poor results in conven- 
tional plotting and triangulation procedures for emitter position fixing. 
The Kalman filtering process, in effect, statistically filters out random 
Signal error or measurement noise, resulting in an optimal estimate of 
ο true signal. 

By computer processing techniques these emitter data may be 
sorted and sequentially filtered to yield optimal estimates of bearing 
angle-of-arrival information. Computerized en as suggested 
within this paper, may then be utilized to compute numerical emitter 
locations in latitude/longitude coordinates to a high degree of accuracy. 
For limited-storage airborne computer applications the Kalman filter - 
ing process and plane triangulation solution should be readily adaptable 
to airborne implementation to yield real-time cockpit display outputs 
with typical estimated emitter position error within three percent of 
the exact value. 

For more precise resolution of emitter location, the spherical 
earth triangulation solution may be considered; or alternately, the 


extended Kalman filtering techniques may be utilized to improve or 
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update plane triangulation solution initial estimates of emitter location. 
These procedures are somewhat more involved, however, and would 


require greater computer memory and processing time. 
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APPENDIX A 
DERIVATION OF SCALAR KALMAN FILTER 
RECURSION EQUATIONS 
The Kalman filter recursion equations from page 11, equations 


(6), (7), and (8), are listed below: 


G() - P&/k- niaoT [noora/x-08097 + Rao ] °? (1A) 
P(k/k) = P(k/k-1) = G(k)H(k)P(k/k-1) (2A) 
P(k+1/k) = O(k+1, k)P(k/k) B(k+1, k)? + Qk) (3A) 


where in this application Ф (к+1, к) = P η , 


H(k) = E o | И 


the P and Q matrices may be considered to have scalar components 


given by 
PIL 5ι2 Q11 Q12 
and respectively, 
τοι P22 (9221 922 | 
G1 


and the G тпаї{гт1х їс а 2х1-тайш 1 хо the form 
G2 


Writing equation (1A) in matrix notation yields 


αι Р ο ρω. E 
= E ο] + R ; (4A) 


G2 BER P22 0 Pal ΕΞ; 0 


and for the observation matrix H(k) = E ο] , the inverse term becomes 
a scalar allowing the gain terms to be computed directly, as is shown 


below: 
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(х1 Pri =l 


= (ри + R | (Б) 
G2 P21 


which results in the scalar gain equations 


P11(k/k-1) | (6A) 


SUD E утку) ἐπα 


P21(k/k-1) 


a u τα oe ° (ТА) 
P11(k/k-1) + R (k) 


From equation (19) which defines the P matrix, it can be seen that since 


8 and д аге statistically independent, where 

~ N 

8 = @ Е 6 , 
then P12 = P21, and the P matrix is symmetric. From equation (15) it 
can be seen that the Q matrix is also symmetric. Equation (7A) then 


becomes 


GNE = ..P12(k/k-1) — 5 ` (8A) 


P11(k/k-1) + R(k) 
To solve for the prediction covariance scalar terms, equation 
(2A) was substituted into (3A) giving 


КОП Р12| [1 ТІРРІІ Р1271 |Сі1 Р11 Р12141 01 (011 012 


: 1. 119] + 
P12 P22] Jo ı |[[pı2 P22] [c2 P12 p22|\irıj |Q012 Q22 
k+1/k k/k-1 k/k-1 


which upon simplifying and writing in scalar form yields equations as 


noted on page 84 of the Kalman filter program. 
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APPENDIX B 


DERIVATION OF SCALAR SMOOTHING EQUATIONS 


The general smoothing equations (11) and (12) are given on page 


13 and are repeated below: 


A A 
X(1/k) = X(1/k-1) + D(1/k)G(k)E(k) (1B) 
where D(1/k) = D(1/k-1)P(k-1/k-1) @(k, k-1)T P(k/k-1)7} (2B) 
A 
and E(k) = Z(k) - H(k)X(k/k-1) . (3B) 


A 
In order to solve for a smoothed value of the state X(1/k), it is first 
necessary to solve equation (2B) sequentially since the state transition 
matrix ® is a function of the variable sampling interval as given by 
: 1 0 
G(kk-1) - (4B) 
T(k) 1 
It is also necessary to obtain the inverse of the error covariance matrix 
P in order to solve equation (2B). On the IBM 360, the matrix inverse 


routine MINV was used with good numerical accuracy, and the results 


were programmed as follows: 


4 [Pu Ῥια]-; [Pmi ріміз 
Р = = ; (5B) 
= P21 P22 PIN21 PIN22 
and letting 
P11 P12 ΠΠ EP Ollie οὐ 
: „Ds πο = 
7 P21  P22 Zl Dy ЕО ar ο» 
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as before. Since the P matrix is symmetric (as shown in Appendix A), 


“lis also symmetric. 


E 
To Seine for the scalar smoothing equations, equation ,(8) on page 

11 was solved for P(k-1/k-1) giving 

P(k-1/k-1) = $0, k-1)7! Е - Qu] o” (S RE (6 B) 


This equation was then substituted into equation (2B) giving 


which written in matrix form yields 


р11 П12 RIED om eo Q11 Q12 pes PIN12| 
D21 D22 ή ο αυ Ὁ 1 012 022) (PIN12 тш) 
1/k bel 


This reduces to 


Dil D12 DII “Delt τι ο ое -QP12 
DI D22 D21 D22 1 -QP21 1.0 - QP22 
la: 1/k-1 


where QP1I = Q11(PIN11) + Q12(PIN12) 

QP12 = Q11(PIN12) + Q12(PIN22) 

QP21 = Q12(PIN11) + Q22(PIN12) 

QP22 = Q12(PIN12) + Q22(PIN22) 
and D(1/0) was initialized as the identity matrix. These are the terms 
used in the Kalman filter program to solve for sequential values of the 


D matrix as seen on page 84. 
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APPENDIX C 


DERIVATION OF PLANE TRIANGULATION 
SOLUTION 


From Figure 1 it can be seen that 


A Aq ~ Ay 


(ап. = шт (1C) 
j 
A Xp 7 Xt 

and tan, = ——T (2С) 
Xp -Yg 


Substituting latitude (L)/longitude (A) coordinates for the X/Y axes of 
the flat earth model results in angular measurement of distance in the 
plane. Since equal angles of latitude and longitude yield equal distances 
of movement only at the equator, a correction factor must be applied to 
longitude angular measure as either pole is approached from the equator. 
At 60° North latitude, for example, one degree of movement 

measures sixty nautical miles while one degree of movement in longi- 
tude only measures thirty nautical miles. This requires that the longi- 
п angular measure be corrected by the cosine oz the local latitude 
coordinate so that the tangent ratio terms will have units of the same 


dimensional size. This results in the equations given below: 


(> - Х.Усові, 

A ` 

“ΠῚ : Ἐάν зс) 
B i 
N n 

tanB, . (Ar A,)cos T (4C) 
Lip - Le 
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These two equations could be solved simultaneously for the two unknowns 
Ly and Àj, except that the unknown Le also appears as a cosine function. 
The method chosen to resolve this difficulty allowed a known midlatitude 
function 


ЕШ 
A a (5C) 


to be used in approximating an intermediate value of emitter target 


latitude TILA where 
; A A 
(Ap - A,)cos(WAV) + L;tanĝ, - L tan 0, 


ШИША = (6С) 
сап; - tan§r 


results from substituting cos(WAV) for cosLy in a simultaneous solution 
of equations (3C) and (4C) for Lyp. It should be noted that this midlati- 
tude approximation is only valid if the aircraft position latitude is within 
a few degrees of the emitter location latitude. This approximation would 
not be valid for HF sky wave radiation transmitted from an emitter at 
great distances from the receiving aircraft. For situations in which the 
approximation is valid, equations (28) and (29) give a corrected solution 
for emitter latitude and longitude based on the intermediate latitude 


function TILA. 
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APPENDIX D 


DERIVATION OF SPHERICAL EARTH 
TRIANGULATION SOLUTION 


Figure 4a describes the geometry of the spherical triangulation 
solution, in which the arc length sides of the spherical triangles must 
be measured from the angle subtended at the center of the sphere. 

Since latitude (L) and longitude (A) lines meet at right angles and 6; 
and 0; are known, then two angles of the spherical triangles are known. 
This allows the side opposite and the side adjacent to the emitter bear- 
ing angle B to be computed. 

Using the Law of Sines and the Angle Law of Cosines for spherical 
triangles, the general solution equations may be derived as follows from 


Figures 4b and 4c: 


sin(A)  sin(B) . sin(B)sin(a) 
. 0 = =—— or sin{A) -—— ы (1D) 

Sine Law: sin(a) sin(b) sin(b) 
Cosine Law: cos(A) = -cos(B)cos(C) + sin(B)sin(C)cos(a). (2D) 
For angle C = 90° equation (2D) reduces to 

cos(A) = sin(B)cos(a). (3D) 
Dividing equation (1D) by (3D) yields 

tan(A) = _sin(A) = tanla) (4D) 


СОБИ 9 sin(b) 
where a and b are angular arcs of the spherical triangles measured at 
the center of the sphere. From Figures 4b and 4c, bis an arc of latitude 


and a is an arc of longitude which must be corrected for the local latitude 
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at which it is measured. Therefore, 


@ 


= (Ат - d)cosLip 


and b=*Lp-L. 


Substituting these values of a and b into equation (5D) results in 


tan [ Ат - A JeosL... | 


tan Ü: 
i 


tan | dep - \ 001, | 
sin(L« - Le) 


and tan 8 f 





(5D) 


(6D) 


(7) 


(8D) 


which give the exact solution for emitter location on a spherical earth 


of constant radius. 
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APPENDIX E 
COMPUTER ITERATION SOLUTION TO 


SPHERICAL TRIANGULATION 
EQUATIONS 


No closed-form solution to the two spherical triangulation 
equations (30) and (3 1) could be found since the unknowns (Lips dp) are 
enmeshed within multiple trigonometric relations. Instead, an itera- 
tive approximation computer solution was developed which allowed 
successive numerical updating of previously computed values. This 
method is only as good as the initial estimate with which it commences, 
so it is best described as a spherical earth correction procedure to the 
plane triangulation solution for emitter location. The procedure is 


concerned with updating values of the following terms 


DTLAI = Lẹ- Lj 
DTLAF = Ln -L, 
DTLOI = An E 
DTLOF * À4-A, 


by sequentially recomputing them from rearranged απ ους of equations 
(30) and (31). The procedure is repeated and successive terms are 
compared until a satisfactory solution is reached. See the computer 
subroutine onthe following page for a further explanation of the iterative 


process. 
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SPHERICAL EARTH SOLUTION SUBROUTINE 
ж ж ak ok ok ake ae ak ote ac aks ЕЗІ ІЗ ЗІ ak ok ake еж: сок ke ok ok ok ok ok ak ake ak ok 


CORRECTIONS TO PLANE TRIANGULATION SOLUTION FOR 
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APPENDIX F 


DERIVATION OF EXTENDED KALMAN FILTER 
SCALAR EQUATIONS 


š The extended Kalman filter recursion equations (42), (43), and 
(44) are similar in form to the basic Kalman filter recursion equations 
on page 11, except that the observation matrix is a non-linear transfor - 


mation matrix M as defined by equation (41) where 


E = | = [омх DMY | (1Е) 
Ro IA AL, i U 


91 
Mk) = 3x 





x 


Taking partial derivatives of equation (24) yields 


(Lip - L.)cosL,, 






DMX . 98 "= —— —  — (2F) 
9 Ary (Lip - L.) + (Ar - №.) cos LT 
д c XV. nec Бб "σος. 

ros 90 _ Т i | T i T - T (35) 


2 
4% (Lip - L4) d: (Ат - Aj σος Lj 


These terms are computed numerically and substituted into equations 
(42) and (43) which are rewritten in matrix form as 


GX ΕἸ ΕΙ ΗΠ Rid 
= (рмх рмү | 


GY IST pm P21 ο hi 


-1 


P12| |DMX 
Ж 


апа 
Pll p121 (1 Ol (eit РОР ОЛ ОО 

š + (SF) 
PIPE ο ρου Ρο 1 Q12 @22 


k+1/k k/k 


a2 





don (43) is then substituted into (5F); however, the fact that the 
o matrix reduces to the identity matrix for the stationary filtered 
states (Lip, dp) considerably simplifies the recursion equations since 
they are no longer dependent on the sampling interval T(k). Since the 
gain terms are computed as a function of the non-linear terms DMX 
and DMY, it is essential to have the correct units on the error term 
E(k). In this case the units of gain are not dimensionless but degrees per 
radian, so the error term must be given in radians to update position 
estimates given in degrees of latitude and longitude. 

These scalar equations may then be written out and programmed 
sequentially as was done in Appendix A for the linear recursion 
equations. See page 92 for a listing of these equations within the exten- 


ded Kalman filter subroutine. 
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APPENDIX G 


MONTE CARLO SIMULATION ΕΠΕῸΣ ВИСИ 
ESTIMATED EMITTER POSITION ERROR IN NAUTICAL MILES 
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AIRCRAFT FLIGHT DATA COMPUTER OUT I T 


LISTING OF BEARING ANGLES-GF-ARRIVAL AND AIRCRAFT NAVIGATION DATA 
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NTAR = 1 


UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED TO NTAR = 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 


NUM 


FREQ 


2872.0 
2872.0 
2872.0 
DET aa 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2812.0 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 


PRF 


320.0 
322.0 
320.0 
318.0 
320.0 
319.0 
321.0 
320.0 
320.0 
320.0 
319.0 
319.0 
320.0 
320.0 
320.0 
320.0 
319.0 
323.0 
320.0 
322.0 


PW 


2.40 
2.50 
2.50 
2.40 
2.60 
2.20 
2.50 
2.60 
2.20 
2.40 
2.60 
2.50 
2.50 
2.40 
2.50 
2.50 
2.30 
2.80 
2.50 
2.50 


HDGD 


135.52739 
135.64830 
135.41739 
135.30740 
135.00026 
135.28560 
135.51651 
135.43913 
135.38472 
131.33047 
121.44293 
120.89346 
121227792 
121.12437 
120.72844 
120.65167 
118.74031 
120.64079 
120.65167 
120.75079 


56 


BRNGD 


283. 69995 
282.10010 
280.80005 
279.90015 
278.49976 
277.30005 
275.89966 
275.40015 
273.79980 
276.00000 
280. 89990 
279.69995 
278.29980 
276.50024 
275.00000 
274.30005 
272.19995 
267.10010 
267.00000 
263.39990 


THETAD 


59. 22 129 
21.174829 
56621729 
55220752 
53.50000 
52+ 58545 
21. 41 602 
50.83911 
49. 18433 
47.33032 
42.342 11 
40.59326 
39.517164 
31.62451 
35.728217 
34.95166 
30.94019 
27.74072 
27.65161 
24.15063 


I 





UNFILTERED EMITTER 


NUM 


3 
20 
24 
30 
35 
43 
44 
50 
61 
65 


FREQ 


2876.0 
2878.0 
2876.0 
2878.0 
2876.0 
το. 0 
2876.0 
2876.0 
ЙЕ 16.0 
BE. 0 


PRF 


306.0 
304.0 
305.0 
306.0 
306.0 
306.0 
306.9 
306.0 
306.0 
305.0 


NTAR = 2 


TARGET DATA INITIALLY ASSOCIATED TO NTAR = 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 


PW 


2.40 
2.40 
2.40 
2.30 
2.30 
2.40 
2.30 
2.50 
2.30 
2.40 


HDGD 


135.48329 
E22, 29930 
120.76172 
121.58500 
129.68495 
119.66331 
Ш 95276851 
118.82796 
121.01439 
119.79510 


ΟἿ 


BRNGD 


286.89990 
291. 9980 
2589 309090 
289.60010 
281.29980 
287.19946 
28509395 
29299915 
27 028972770 
ο ο ο το 


THETAB 


62.38306 
54.09912 
19.161682 
912165005 
4+7. 98462 
4+6. 86255 
44, 978093 
41. 82764 
38.91 406 
31.894718 


2 





NIAR = 3 


UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED TO ΙΤ... 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 


NUM FREQ PRF PW HDGD BRNGD THETAD 
5 2846.0 370.0 2:30 Pes boos 300.99976 16. 31836 
18 2846.0 310.0 2.0 134.45079 299.49976 13.95044: 
29 2845.0 370.0 2.40 121.57414 308.69971 10.217441 
34 2846.0 Ste) 2.40 119.94865 309.39990 69. 34839 
41 2846.0 370.0 2.40 120.66260 306399976 61.66235 


58 





NTAR = 4 


UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED TO NTAR = 


EMITTER BEARING LINES GF SIMILAR FREQUENCY AND PRF OR PRE 


МОМ 


13 
16 
23 
25 
27 
36 
&6 
49 
58 


64 


FREQ 


5556.0 
5554.0 
5556.0 
5556.0 
5556.0 
5554.0 
5554.0 
5554.0 
5554.0 
5556.0 
5556.0 
5556.0 


PRE 


320.0 
640.0 
320.0 
640.0 
320.9 
640.0 
640.0 
640.0 
320.9 
640.0 
640.0 
320.0 


PW 


1.00 
1.10 
1.90 
1.00 
1.10 
1.10 
1.10 
0.90 
1.00 
1.00 
1.00 
1.10 


59 


HDGD 


135.06616 
135.49416 
135.30740 
120.76172 
120.87172 
121.36617 
120.70670 
118. 83884 
118.72943 
120221225 
121.15704 
119.75102 


BRNGD 


314. 80005 
312.29956 
311.60010 
276.30005 
275.19995 
275.09985 
272.19995 
312.00000 
311.09985 
264.00000 
251. 09985 
251.00000 


THETAD 


89.86621 
81. 19510 
86. 90 (47 
37, 06177 
36.071553 
36.46582 
32.90649 
{0.83862 
69.82910 
24.21216 
18. 22684 
16. 75098 


А 


MULTIPEE 





NTAR = 5 


UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED TO NTAR = 5 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 


NUM 


11 
17 
21 
26 
38 
45 
a 
55 
60 


FREQ 


5506.0 
5504.0 
5504.0 
5504.0 
5504.0 
2506.0 
5504.0 
5504.0 
5504.0 


PRE 


640.0 
640.0 
640.0 
640.0 
640.0 
640.9 
642.0 
640.0 
642.0 


PW 


0.60 
0.60 
0.60 
0.60 
0.60 
0.T0 
0.60 
0.60 
0.60 


HDGD 


135.61504 
135.16473 
120.58583 
120.98116 
121.17940 
1178799298 
119.64192 
120.88258 
120.77260 


60 


BRNGD 


312.19980 
3211. 29980 
2194609971 
dS ROO 
315.00000 
312.79980 
309. 00000 
306. 10034 
3201 59966 


THETAD 


88.41479 
86.464236 
36.28540 
35.18066 
76.17920 
71.792772 
68.64087 
66.98291 
62.17212 





NTAR = é 


61 


UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED TO NTAR = 6 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRE OR ΡΑΕ ποτ. NS 
NUM FREQ PRE PW HD GD BRNGD THETAD 

ШІ 28/6.0 23970 2.20 135. 56006 21 {ο 30005 925 860 

22 2876.0 359.0 2.20 120.73991 283. 30005 44. 03979 

47 2876.0 359.0 2.40 118. 71851 216.50024 35.21515 

52 2874.0 32920 2.40 120.40988 ¿(1.899990 32902501 

56 2876.0 299.0 2.20 Dea | 266280005 21.92432 

62 2874.0 359.0 2.40 121.15704 2986. 59966 59.55664 





NTAR = 7 


UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED ID DT A P P DE 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 
NUM FREQ РКҒ РИ HDGD BRNGD THETAD 


31 9268.0 500.0 0.20 120.76172 287260010 48.36182 
39 9268.0 499.0 0.20 120.91527 283559966 44.31470 
57 9279.0 500.0 0.30 121.22295 2(6.59985 212582212 


62 





DATA 
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NT 


FILTERED AND SMOOTHED EMITTER 


NUM 


NO PN om 


9 
10 
14 
15 
19 
28 
32 
33 
37 
40 
42 
48 
53 
54 
59 


SMOOTHED INITIAL BEARING ANGLE 
FILTERED FINAL BEARING ANGLE 


FREQ 


2872.0 
237220 
2872.0 
2872.0 
2872.0 
87250 
2872.0 
2872.0 
2872.0 
2872.0 
812.20 
2872.0 
2872.0 
2872.0 
2872.0 
2872.0 
2572.0 
2872.0 
2572.0 
2872.0 


PRE 


32070 
322.0 
320.0 
318.0 
320.0 
319.0 
321.0 
320.0 
320.0 
320.0 
219.0 
31930 
220.0 
320.0 
320.0 
320.0 
319.0 
σον ο 
320.0 
52230 


PW 


2.40 
2.50 
2.50 
22.40 
2.60 
2.20 
2.50 
2.60 
2.20 
2.40 
2.60 
2.50 
2.50 
2.40 
2.50 
2.50 
2.30 
2.80 
2.50 
2.50 


AR = 1 


TARGET DATA CORRELATED TO NTAR = 


THETAD URTO 
νο ο ο NZ] 
51.14829 51.14831 
56. 217129 26612552 
55220752 95215414 
53.50000 532292235 
52.586545 52.52164 
51.41602 51. (4019 
ΠῚ 50.47720 
49.18433 49.66541 
417233032 47281853 
а2.54277 42. 15356 
&0.59326 40.52049 
39.57764 39. 117419 
37.6245 1 Ste 25296 
33. 12827 35.84486 
34.95166 35.12247 
30.94019 31207832 
27. 14072 28. 14043 
21.65161 21. 15360 
24215063 24.710815 

= 59.56096 
= 24.10815 


PLANE TRIANGULATION SOLUTION OF EMITTER LOCATION 
ЕПІТТЕК EAT Tove 
EMITTER LONGITUDE 


5 


66 


3.98666N 


=-119.60387W 


THT D1 


5922223 
Peco veo 
59. 01391 
59.04889 
59.18648 
59.16264 
59.24341 
59.14085 
59.24218 
59.36894 
59. 51895 
59.50125 
59. 53 006 
59.51884 
59.52922 
59.53961 
59.54620 
59.55690 
59.55646 
59.56096 


1 





NT 


AR = 2 


=I TERED AND SMOOTHED EMITTER TARGET DATA CORRELATED ΤῸ ΠΠ 


NUM 


3 
20 
30 
35 
43 
44 
50 
61 
65 


SMOOTHED INITIAL BEARING ANGLE 
FILTERED FINAL BEARING ANGLE 


FREQ 


2876.0 
2878.0 
2878.0 
2876.0 
2876.0 
2816.0 
2816.0 
2816.0 
2876.0 


PRF 


306.0 
304.0 
306.0 
306.0 
206.0 
306.0 
306.0 
306.0 
306.0 


PW 


2.40 
2.40 
2.30 
2.30 
2.40 
2.30 
2.50 
2.30 
2.40 


THETAD THTD 
62.38206 62.28306 
54.09912 54.09912 
51.128555 51.031205 
41.984562 48. 03497 
46.86255 46.09224 
44.971803 45.39227 
41.82 [6% 43.05530 
38.91406 оо е у Б 
37289478 31.80147 

= 62.21768 
= 37.80147 


PLANE TRIANGULATION SOLUTION OF EMITTER LOCATION 
EMITTER LATITUDE 
EMITTER LONGITUDE 


3 


67 


4.12234М 


=l 19700333W 


THTD1 


62. 38306 
62.38306 
62: 29199 
62. 32124 
61.92012 
62. 0042 f 
62.22186 
62.22676 
62.21768 
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NTAR = 3 


=ILTERED AND SMOOTHED EMITTER TARGET DATA CORRELATED TO NTAR = 
NUM FREQ PRE PW THETAD THTD THTD1 

ИЕЛІ 270200 220 76.31836 76.318236 το τα, 
JI 2846.0 370.0 2.20 73.95044 — 73.95042 76.31836 
Bo 2846.0 370.0 2.40 70.27441 70.56093 76.49266 
КО 2846.0 370.0 2.40 69.34839 69.21106 76.45313 
mm 2846.0 370.0 2.40 67 «66235 67.62057 16.44316 
SMOOTHED INITIAL BEARING ANGLE = 76.44316 

FILTERED FINAL BEARING ANGLE = 67.62057 


PLANE TRIANGULATION SOLUTION OF EMITTER LOCATION 
EMITTER LATITUDE = 33.81570N 
EMITTER LONGITUDE =-118.156428 
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FILTERED AND SMOOTHED 


NUM 


8 
13 
16 
46 
49 


SMOOTHED INITIAL BEARING ANGLE 
FILTERED FINAL BEARING ANGLE 


FREQ 


5556.0 
5554.0 
5556.0 
5554. 0 
5554.0 


PRF 


320.0 
6&0.0 
320.0 
640.0 
320.0 


NTAR = 4 


EMITTER TARGET DATA CORRELATED TO NTAR = 


PW 


1.00 
1.10 
1.00 
0.90 
1.00 


THETAD THTD 
89.806621 89. 86621 
87.79370 87.725370 
86.590747 86. 86319 
70.83862 70.91100 
65. 82910 69. 16620 

= 90.18146 
= 69.16620 


PLANE TRIANGULATION SOLUTION OF EMITTER LOCATION 
EMITTER LATITUDE 
EMITTER LONGITUDE 


33.25166N 


69 


=-119.49638W 


ΤΗΤΟΙ 


89.86621 
89.86620 
89. 84418 
90.19025 
90.18146 
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ЕСЕТЕРЕО АМО SMOOTHED EMITTER TARCET DATA 


NTAR = 5 


CORRELATED TO NTAR = 


NU M FREQ PRF PW THETAD THTD THTD1 
11 5506.0 640.0 0.60 88.41419 88.41479 88.41479 
17 5504.0 640.0 0.60 86.46436 86.46436 88.41478 
38 5504.0 640.0 0.60 πο ιτ 16.17416 88. 28852 
45 5506.0 640.0 0.70 (1. 75212 (2.21440 88.49161 
51 5504.0 642.0 0.60 68.64087 69. 23541 88. 57054 
55 5504.0 640.0 0.60 66.98291 66.89281 88.56223 
60 5504.0 642.0 0.60 62.11212 63.35942 88.63921 
SMOOTHED INITIAL BEARING ANGLE = 88.63921 

FILTERED FINAL BEARING ANGLE = 63.35942 


PLANE TRIANGULATION SOLUTION OF EMITTER LOCATION 
EMITTER LATITUDE = 33.25288N 
EMITTER LONGITUDE =-119.45715W 
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NTAR = 6 


FILTERED AND SMOOTHED EMITTER TARGET DATA CORRELATED TO NTAR = 


NUM 


12 
22 
47 
52 
56 


SMOOTHED INITIAL BEARING ANGLE 
FILTERED FINAL BEARING ANGLE 


PLANE TRIANGULATION SOLUTION 


FREQ 


2876.0 
2876.0 
2876.0 
2874.0 
2876.0 


PRE 


359.0 
359.0 
359.0 
359.0 
359.0 


PW 


2.20 
2.20 
2.40 
2.40 
2.20 


EMITTER LATITUDE 
EMITTER LONGITUDE 


THETAD THTD 
52.86011 52.86011 
4.03979 t4. 03979 
293321875 34.57146 
32.309517 32.326 140 
24.92432 28.921722 

= 52.26216 
= 28.921722 


OF EMITTER LOCAT ION 


33.99094N 
A бү 


(1 


THTD1 


52.86011 
52.86009 
52204729 
52 05292 
52.26216 
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NTAR = 7 


FILTERED AND SMOOTHED EMITTER TARGET DATA CORRELATED TO NTAR = 7 


NUM FREQ PRE PW THET AD THT D THTD1 
ЕШ 9268.0 500.0 0.20 48.36182 48.36182 48.36182 
29 9268.0 499.0 0.20 44.31470 44.31470 48.36180 
ЖИ 92710,0 500.0 0.30 31.82275 ο πι 95 4+8. 22450 
SMOOTHED INITIAL BEARING ANGLE = 48. 22450 


FILTERED FINAL BEARING ANGLE Эт. 23 


PLANE TRIANGULATION SOLUTION OF EMITTER LOCATION 
EMITTER LATITUDE = 34.14626N 
EMITTER LONGITUDE =-119.10936W 
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MONTE CARLO SIMULATION COMPUTER OUTPUT 


LISTING OF BEARING ANGLES-OF-ARRIVAL AND AIRCRAFT NAVIGATION DATA 
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-— 
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DATA INITIALLY ASSOCIATED TO NTAR 


EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 


UNFILTERED EMITTER TARGET 


BRNGD THETAD 


PRE PH HDGD 


FREQ 


NUM 
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UNFILTERED EMITTER TARGET DATA INITIALLY ASSOCIATED TO NTAR 
EMITTER BEARING LINES OF SIMILAR FREQUENCY AND PRF OR PRF MULTIPLE 


BRNGD THETAD 


PRE PW HDGD 


FREQ 


NUM 
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